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One of the most important problems in hybrid systems is the reachabil- 
ity problem. The reachability problem has been shown to be undecidable 
even for a subclass of linear hybrid systems. In view of this, the main focus 
in the area of hybrid systems has been to find effective semi-decision pro- 
cedures for this problem. Such an algorithmic approach involves finding 
methods of computation and representation of reach sets of the continu- 
ous variables within a discrete state of a hybrid system. In this paper, 
after presenting a brief introduction to hybrid systems and reachability 
problem, we propose a computational method for obtaining the reach sets 
of continuous variables in a hybrid system. In addition to this, we also de- 
scribe a new algorithm to over-approximate with polyhedra the reach sets 
of the continuous variables with linear dynamics and polyhedral initial 
set. We illustrate these algorithms with typical interesting examples. 

1 Introduction 

Hybrid systems combine discrete and continuous dynamics. The dynamics of the 
continuous variables within a discrete state are specified by differential equations 
or differential inclusions. An important problem in the analysis and sysnthesis 
of such systems is the so called reachability problem, which asks, for two sets of 
configurations of a given hybrid system, say X\ and A2, where a configuration 
consists of discrete and continuous components, whether or not there is a hybrid 
trajectory with the initial configuration in X\ and the final configuration in X^. 
A hybrid trajectory may be described as a trajectory of configurations consist- 
ing of discrete state jumps and smooth arcs, where each arc evolves according to 
the continuous dynamics of a discrete state, with the starting and end points of 
each arc satisfying the jump conditions of discrete transitions. A more precise 
definition is presented in Sec. 2, which provides a concise introduction to hybrid 
systems and the reachability problem. 
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The reachability problem is undecidable for certain classes of linear hybrid 
systems (i.e., hybrid systems having linear trajectories within each discrete state 
and linear or constant reset maps, also called constant slope hybrid systems) 
[0, |2) , although in some cases decidability results have been obtained || || (see 
also Q). Therefore, for a general hybrid system, a reasonable alternative a 
pears to be to find semi-decision procedures for the reachability problem [^, ^, 
A computational approach to this problem also requires finding the set of states 
reached by the continuous variables evolving according to the dynamics of a 
discrete state. In this paper, we consider the problem of computing and repre- 
senting the reach sets of the continuous variables within a discrete state when 
the dynamics of the continuous variables are specified by differential equations 
with initial conditions belonging to a specified initial set. 

In this context, various methods have been proposed in the literature for 
finding reach sets of continuous variables || ^, 10 . In Sec. 3, we describe a 



method for computing the reach sets based on the idea that a subset of the 
boundary of the initial set may be found, such that it is sufficient to compute 
the solutions with the initial points lying in that set. This method is similar 
to that in || (and also to some extent to that described in Q). In particular, 
we present a schematic algorithm, which is somewhat more general in its scope 
than that in Q and simpler compared to that in 0. 

Besides these algorithms for reach set computation, an equally important 
issue the representation of the reach sets for manipulating the sets efficiently. 
This requires representation of the reach sets in terms of more convenient sets 
such as, for example polyhedral or subalgebraic sets, that are simple to represent 
and easy to handle for practical purposes. However, since the representing class 
of sets may not contain a member that exactly equals the reach set, we may have 
to find approximations to the reach sets by those that belong to the representing 
class. (See, for example, |l^, |l2j for approximation with polyhedra, and [jllj for 
approximation using ellipsoids). Typically, over-approximations may be used 
for verifying whether a safety requirement may potentially be violated by any 
of the behaviours starting from a given initial set, while under-approximations 
are needed for characterizing a set of states from which a desirable property is 
always achievable. We describe in Sec. 4 a method for over-approximating the 
reach sets by polyhedra when the dynamics are specified by linear differential 
equations and the initial set is a polyhedron. These algorithms are illustrated 
with some simple examples in Sec. 5, while Sec. 6 concludes the paper. 

2 A Brief Survey of Hybrid Systems 

In this section, we present a brief introduction to hybrid systems, and provide 
motivation for the remaining sections. 

2.1 Preliminary Definitions 

We begin with a somewhat detailed but general definition of hybrid systems. 
Definition 1 A hybrid system is a tuple H = (Q, X 7 X, G, E, Init, f), where 
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1. Q is a finite set of discrete states (also called locations,). 

2. X = R™, n > 1, is the set of continuous states, where R denotes the set 
of real numbers. We demote the continuous state variable by x. 

3. E is a finite set of discrete events or environment actions. Some events 
in E are controllable, while others not. Hence, it is convenient to assume 
that E = E c IJ E^, where E c is the set of controllable events, and E^ is the 
set of uncontrollable or disturbance events. 

4- G C Q x X is a set of state invariance conditions. When the system is in 
state q, the continuous variables belong to G(q) = {x G X : (q,x) € G}. 

5. E C Q x V{W l ) x E x {R™ -> R"} x Q is the set of transition edges. An 
edge e G E, where e = (q e , X e , er e , r e , q' e ), is interpreted as follows: 

• // the continuous state is in X e and the event o~ e occurs, then tran- 
sition edge e is enabled in state q e . Thus X e is the set of switching 
points of the continuous variables from state q e to q' e . The set X e is 
specified by a predicate, and is also called a guard. In case there are 
many transition edges simultaneously enabled, then the system may 
select one of the edges nondeterministically . 

• If a transition edge e is selected by the system, then the continuous 
state is reset using the function r e when the system enters state q' e . 
The reset values obey the state invariance condition, r e (X e ) C G{q' e ). 

6. Init C G is a set of initial conditions. 

7. / is an n-dimensional vector field with real-valued components governing 
the dynamics of the continuous state x. The domain of definition of the 
function f is a set T>, where T> — Q x X (if there are no continuous control 
variables) orV = QxXxU (if U C R m ; m > 1, is the range of the 
continuous control variable, which is a function u : [0, T q ) —> U, where 
T q > is as large as may be necessary, depending on the discrete state 
q € Q). When in state q, the continuous variables evolve according to the 
dynamical law 

§-'<*"» 

or according to the dynamical law 

f t =/(<?, *,u) 

depending on the presence or absence of the continuous control variables!^ 

The initial conditions of the continuous evolution are specified either by 
the reset maps when the system takes a discrete transition and reaches the 
state q, or by the (nondeterministic) initial conditions of the system given 
as the set Init. 



1 Another possibility is to specify the system dynamics of the continuous variables in terms 
of differential inclusions. In this case, f(q, x) is a set, and the continuous variables evolve 
according to 

£ 6 /<«'*>■ 
dt 

But we do not consider this case here. 
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//, in addition to the set of initial values, Init, the set of final or accepting 
values in Q x X , denoted by Final, are also specified, then the hybrid system is 
called a hybrid automaton. 

Referring to the system dynamics in the above definition, for our purposes, 
we assume that there are no continuous control variables. Therefore, in the 
sequel, we assume that the dynamics of the continuous variables are specified 
in the form 

§-'<•■"» 

without any continuous control variables. 

Definition 2 Let H = (Q, X, E, G, E, Init, f) be a hybrid system. Each point 
(<7,x) G Q x X is called a configuration of the hybrid system H. 

In Definitions |^ - |^, we introduce certain terminology, starting with the no- 
tion of a step of a hybrid system, leading upto the notions of predecessors and 
successors of a configuration (q, x). But before that, we mention that predeces- 
sors and successors may be defined irrespective of initial conditions, whereas for 
defining the executions of a hybrid system, we require the specification of the 
initial set Init, as will be seen later (refer to Definitions ||, || and ^) . 

Definition 3 Let H = (Q, X,T,,G, E, Init, f) be a hybrid system, and let 
(q, x) and (g',x') be two configurations of H . Then, the pair of configurations 
((q, x), (q' , x')) is called 

1. a time-step, if q — q' and, for some t > 0, there is a function y : [0,t] — > 
X, satisfying y(s) = f(q,y(s)), y(s) G G(q), for s G [0,t], y(0) = x and 

y(*) = x'. 

2. an edge-step if there is an edge e — (q e , X e , a e , r e , q' e ) £ E, such that 
q = q e , q' = q' e , x G X e and x' = r e (x). 

3. a cr-step, where a G E, if there is an edge I = (qi, Xe,o~e,re, q' e ) € E with 
a = ag, such that q — qe, q' = q'e, x € Xe and x' = r£(x). 

A step of the hybrid system H is a pair of configurations ((q,x), (q' such 
that ((q, x), {q' , x')) is either a time-step for some t > 0, or an edge-step for 
some e S E, or a a -step for some a G E. 

Of course, every edge-step is a cr-step for some a G E, and conversely, every 
(T-step is an edge-step for some edge e G E. Hence these definitions may seem 
redundant. However, the distinction between the two types of transitions would 
become obvious if there is a discrete state controller, defined as a function C : 
Q x X — > E c , that triggers discrete controllable events to enable a particular 
transition depending on the present configuration (see for instance |l3| for the 
case of timed automata). But in this work, we shall not have occasion to discuss 
about discrete state controllers. 

Definition 4 (Trajectories and Executions of a Hybrid System) A hy- 
brid trajectory or simply a trajectory of H is a sequence of configurations 
(<7i,xi), (q 2 , x 2 ), (<?3,x 3 ), ■ • ., where for each i > 1, ((q h Xj), (q i+ i, x i+1 )) is a 
step. A trajectory is 
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1. finite, if the number of steps is finite; 

2. an execution, if (qi,Xi) € Init, where (gi,xi) is the initial configuration; 
and 

3. a finite execution, if it is both. 

We now define, for each a G E, the set-valued successor function, Post a : 
Q x X — > 2*2 x 2 X . As will be briefly mentioned later (after Definition ^) , our 
definitions lead in a natural way to extract a labeled transition system (see GJ) 
from a hybrid system. 

Definition 5 (cr-Successors of a Configuration) 

Let H = (Q, X, E, G, E, Init, f) be a hybrid system, a S £ and (g, x) S Q x X . 
Then, we define Posi CT (g,x) to be the union of the two sets Si,S% Q Q x X, 
w/iere 

1. Si — {(g',x') : ((g,x), (g',x')) is a time-step for some t > 0}; and 

£ = {(g",x") : ((g',x'),(g",x")) «s a a-siep /or some (g',x') € Si }. 

Further, j/(g s ,x s ) € Post a (q,x.), then (g s ,x s ) is called a cr-successor o/(g,x). 

The set Si, as in the first part of this definition, is of main interest to us in the 
later sections of this paper. Specifically, we will deal with the set Reachf (X ), 
where g € Q and Xq C G(q), defined as 

Reach q (Xq) — {x : ((g, xq), (q, x)) is a time-step for some t > and Xo € Xo } . 

(1) 

In the subsequent sections, we will be concerned more with this and a related set. 
But for now, we shall proceed with our discussion with the following definition: 

Definition 6 (cr-Predecessors of a Configuration) 

Let H = (Q, X, E, G, E, Init, f) be a hybrid system, a G E and (g, x) G Q x X . 
Then, we define Pre CT (g,x) to be the union of the two sets Pi,P^ Q Q x X, 
where 

1. Pi = {(g',x') : ((g',x'), (g,x)) is a time-step for some t > 0}; and 

2. P 2 CT = {(<?", x") : ((g",x"), (g',x')) is a a-step for some (g',x') e P l }. 

Further, if(q p ,Xp) £ Pre a (q,~x), then (g p ,x p ) is called a cr-predecessor o/(g,x). 

Note that the notion of tr-successor generalizes the notion of tr-step by includ- 
ing time-steps. The set valued function Post a (q,x.) defines, in a natural way, 
a transition relation, -^-> C (Q x X) x (Q x X), as follows: (q, x) (g',x'), 
if (g',x') £ Post a (q,x). The transition relation is also called a -transition 
relation. Similarly, the set valued function Pre CT (g,x) defines another tran- 
sition relation, C (Q x X) x (Q x X), as follows: (g, x) <-^— (qf, x'), if 
(g,x) e Pr e<T (g',x'). 

For two configurations (gi,xi) and (g 2 ,x 2 ), if ((gi,xi), (g 2 ,x 2 )) is a er-step, 
then (g 2 ,x 2 ) S Posi(gi,Xi) and (gi,xi) € Pre(g 2 ,x 2 ); hence, in this case, we 
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have (gi , Xi ) (g 2 , X2 ) and (qi , xi ) (g 2 , X2 ) . This may mislead us to get the 
false impression that — —> and are the same; and to avoid any such possible 
confusion, we emphasize that, in general, it is not true that [q%, xi) — > (g>2, X2) 
implies (qi, xi) <-^— (52, X2); and the same with the converse statement. Hence, 
the two relations — and «-^— are not the same. 

The definitions of Post a and Pre a can be extended straightforwardly to 
subsets ofQxX, as follows: for 5 C Q x X, define 

Post a (S) = (J Post a (q,x), and Pre ff (5) = |J Pre a {q,x). 
(?,x)es (g,x)es 

Finally, define 

Pos<(5) = |J Post a (S), and Pre(S) = (J Pre a {S). 

We sometimes refer to Post and Pre as the 1-step transition functions^]. More 
generally, the k-step transition functions Post k and Pre k are defined inductively 
as follows: For S C Q x X, 

PosfiS) = Post(S), and for A; > 2, Post k (S) = Post(Post k - 1 (S)). 

Similarly, 

Pre 1 (S) = Pre(S), and for k > 2, Pre k (S) = PreiPre^ 1 (S)). 
Finally, define Post* and Pre* as 

Post*(S) = (J Post k (S), and Pre*(S) = |J Pre k (S). 

k>l k>l 

2.2 Reachability Problem for Hybrid Systems 

In the notation just discussed, the reachability problem for a hybrid system 
H = (Q, X, S, G, E, Init, /), where Init is not necessarily specified in advance, 
may be posed as follows: 

ReachProbleml: For two subsets Si and S% of Q x X, is there a finite 
trajectory, (q 1: xi), (q 2 , x 2 ),..., (<7jv, xjv), for some 
N >1, such that (gi,xi) € Si and (<7jv,Xjv) e5 2 ? 

This may also be rephrased as 

ReachProblem2: For two subsets Si and S2 of Q x X, 
whether Post* (Si) f| S 2 ^ ? 
equivalently: whether Si f] Pre* (S2) 7^ ? 

If the initial set Init is specified in advance, then, with Si = Init, the reacha- 
bility problem ReachProbleml becomes 



2 This should not be confused with the notion of step as denned in Definition 
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ReachProblem3: For a subset S of Q x X, 

is there a finite execution of H 
with final configuration in S ? 



There is a counterpart to the reachability problem, called the avoidance 
problem, which may be posed as follows: 

AvoidProblem: For two subsets Si and S2 of Q x X , 
whether Post* (Si) f| £2 = ? 
equivalently: whether Si f] Pre*(S 2 ) = ? 



In the sequel, we restrict our attention to ReachProbleml or ReachProb- 

lem2. We observe that the answers to these questions depend, in general, not 
only on the hybrid system H, but also on the class C of subsets of Q x X that 
are under consideration. Informally, the class C is required to 

1. include a specified class of sets S, where S may consist of the initial 
set, sets defined by state invariance conditions, those defined by guard 
conditions, and domains and ranges of reset maps of edges (and also final 
sets, if specified), 

2. be closed under the boolean set operations of union and complimentation, 
and under the functions Post a and Pre a , and 

3. have an effective decision procedure for answering questions such as, for 
two sets Si, S2 € C, whether Si = S2 or not. 

A formal presentation of these notions is beyond the scope of this work, 
although a brief discussion may be found in Appendix A. For more details, 
the interested reader may refer to the references on Model Theory, such as 
[fl9[ [20| , j2j| , p2|| . But, for our purposes, we will be content with the following 
(somewhat informal) definition: 

Definition 7 An algorithm for the reachability problem is said to be 

1. a decision procedure, if the algorithm stops after a finite number of steps 
with the correct answer, where the answer can be either yes or no. 

2. a semi-decision procedure, if the algorithm 

(a) never stops with an incorrect answer, and 

(b) always stops after a finite number of steps with the correct answer, 
whenever the answer is yes. 

A hybrid system H is decidable, if there is a decision procedure for the reacha- 
bility problem for H . 

Obviously, for a hybrid system H , if there are semi-decision procedures for 
both reachability problem and avoidance problem, then H is decidablc. 

We describe below a schematic semi-decision procedure for the reachability 
problem of a hybrid system: 
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Semi-Decision Procedure for ReachProblem2 



Input : Sets Si and S2 

Output : "yes" , if there is a trajectory from Si to 52 



Initialization S := Si 

while Sf]S 2 = do 

S := Post{S) 
end while 



return "yes" 

In the later part of this section, we shall be mainly concerned with the 
computational aspects of the Posi-operator appearing in the while-loop of the 
above schematic. 



2.3 Computation of the Post-operator 

We describe here a computational approach to finding Post(S), S C Q x I 
as in the semi-decision procedure for ReachProblem2 discussed in Sec. 2.2 



Recall that, within a discrete state q £ Q, the continuous dynamics of x are 
specified by 

§=><«.*> < 2 > 

with the initial conditions x(0) = Xo £ Xo C X = R n . We begin our discussion 
with the following definition. 

Definition 8 Let R + denote the set of nonnegative real numbers. Then, for 
each discrete state q £ Q, the flow associated with equation is a function 
4> q : R n x R+ -> W l , defined as <j) q (x,t) = ~/ x {t), t > and x £ R", where the 
function 7 X : M + — > R™ satisfies 

—jj— = /W.7x(*)), 
with the initial condition 7 X (0) = x. 

Recall the set of reachable phases Reach^ as defined in (|l|) . In the following 
definition, the same is defined in terms of the flow function cj) q : 

Definition 9 Let q £ Q and Xq C R™. Then, the set of reachable continuous 
phases in the state q, is the set Reach^ (Xq) defined as 

Reach^(Xo) = {x £ R™ : 3x £ X 3t > x = g (x o , t) and 

Vs < s < t <^(x ,s) £ G(q)} (3) 

We now define another operator, called the projection operator, as follows: 

Definition 10 Let S C Q x X . Then, the projection of S on X, is the set 
ir x (S) = {xeX : 0,x) £ S}. 



In terms of the operators ir x and Reach^ , the operator Post may be rewritten 
as 

Post(S)= (J {q' e } x [r e o Reach q o n x (S)], (4) 

e—(q e ,X e ,a e ,r e ,q' e )£E, q—q e 
X e fl [ReachG o tt x (S)] ^0 

where, for two operators g and h, g o h is the composition of the two operators. 
With regard to computational issues, the set operators tt x and r e may not pose 
difficulties (provided r e is suitably specified). Hence the problem reduces to 
that of computation of Reach^, and effectively representing the resulting set. 
Referring to the definition of Reach^ as in (^), the significant and challenging 
task is the elimination of quantifiers, whenever possible. But, as shown in p3fl , 
not all theories of the real number system may admit quantifier elimination (see 
also the discussion presented in and in Appendix A). This fact provides a 
motivation for the study of alternative methods (without requiring quantifier 
elimination method) for computating the set Reach^(Xo). 

Specifically, in the remaining part of the work, we shall be concerned with 
computation of the operator Reach 1 ^ and with representation of the resulting 
sets. Below is a concise description of the questions that we shall be interested 
in and of the organization of the rest of the paper: 

f. How to effectively compute the set operator Reach q , at least when the 
global invariance condition is not imposed but a time bound is specified? 
In this case, we are interested in the set Reach q (XQ, [0, r]), r > 0, Xq C 
R", defined as 

Reach q (X Q , [0, r]) = {xel": 3x G X 3t < t < r and x = 9 (x o , t)} . 

In Sec. 3, we present a schematic algorithm for computing Reach q (X , [0, r]) 
This algorithm is based on a generalization of the method described in || . 

2. How to compute - either exactly or approximately - the set operator 
Reachg when the global invariance condition is imposed (without time 
bound). In this case, we have to deal with the set operator Reach q as 
defined in (^) , and obviously, it would be best to find an algorithm based 
on quantifier elimination. But, as mentioned earlier, we do not assume 
that quantifier elimination method is feasible, hence we may have to find 
an algorithm that computes an approximation to Reach^ . In Sec. 3, we 
indicate how to extend the algorithm for computing Reach q (Xo, [0, r]), 
as mentioned above, to compute an under-approximation or an over- 
approximation (depending on which one is preferred) to the set Reach^ (Xq) . 

3. Finally, how to represent the reach sets obtained by the operators Reach q 
and Reach q , i.e. the sets Reach q (X , [0, r]) or Reach^(X ) as the case 
may be, such that boolean set operations such as union and complementa- 
tion can be performed efficiently. This requires represention of the reach 
sets in terms of sets that are simple and easy to handle, such as poly- 
hedral sets and subalgebraic sets. In Jl2| ], an algorithm algorithm for 
over-approximating the reach sets using polyhedral sets is presented. In 
Sec. 4, we describe another algorithm for over-approximating the reach 
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sets with polyhedral sets. The algorithm presented in this paper differs 
from that of pl| , but if these two methods - that described in Sec. 4 and 
that of Q - are used together, then better results of over-approximation 
may be obtained. 

The remaining part of the paper is orgainized as follows. In Sec. 3 and Sec. 
4, we restrict our attention to these issues. We illustrate these algorithms with 
simple examples in Sec. 5, while Sec. 6 summarizes the paper. 

3 Computational Approaches for Finding Reach 
Sets 

3.1 Preliminary Discussion on Reach Sets 

Let H = (Q, X,Y,,G, E, Init, /) be a hybrid system. For q S Q, consider the 
equation 

f = /(,,*) (5) 

with the initial conditions x(0) = xo £ Xo C X = R n . It is customary to assume 
that f(q, •) is Lipschitz continuous in the second variable, in order to ensure 
existence and uniquness of solutions to the system of differential equations in 
(||). Further, we assume that the initial set Xo is closed. Let, as before, the flow 
(refer to Defintion ||) associated with (||) be <j) q (x,t) and let Reach q (X ,t) — 
{(/) g (xo, t) : xo £ Xq} be the set of phases reached at time t > in state q with 
initial conditions in Xq- Further, let 

Reach q (Xo, [0, r]) = Reach q (Xo,t), r > 0, 

0<t<T 

Reach q (Xo, [0, oo)) = Reach q (Xo, [0, r]). 

r>0 

We shall now define, for a set X q , such that Xo C X q C M™, another set 
Reach' q (X ,X q7 [0, oo)). For this purpose, we first define, for each x £ X a , a 
number r x as follows: 

r x = inf{t >0 :0 ? (x,t) £ X q }, (6) 

where we assume that if r x = oo, then </> g (x, oo) denotes the set of all w-limit 
points (see, for instance, of <fi q (x.,t), for t > 0, and 

{0 9 (x, t) : < t < oo} = {(f) q (x, t) : < t < oo} |J </>(x, oo). 

Further, let 6 g (x, X q ) C W 1 be defined as follows: 

{{0 g (x, t) :0<t< t x }, if r x < oo and </>(x, r x ) E X q , 
{4> q (x 7 1) : < t < oo}, if t x = oo and 0(x, oo) C X g , 
{</> g (x, t) : < t < t x }, otherwise. 

Then the set Reach' q (Xq, X q , [0, oo)) is defined as 

Reach' q (X , X g , [0,oo)) = |J e,(x,X g ). 

xex 
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Now, by taking X q — G(q) in the above definition, we have 

Reach^ (Xq) = Reach' q (Xq, X q , [0, oo)) , 

where Reach^ (Xq) is defined as in (|^), reproduced below for convenience: 

Reachf(X ) = {x G R" : 3x e X 3t > x = g (x o , t) and 
Vs < s < t g (x o , s) e G(q)} 

Throughout the rest of the discussion, we fix a state q, and consider the 
problem of computing the sets Reach q (XQ, [0, r]) and Reach' q (Xq, X q , [0, oo)). 

3.2 A Computational Method: Generalized Face-Lifting 
Algorithm 

When Xq is suitably specified (such as, for example, a rectangle), this prob- 
lem may be solved by finding the solutions to (^) with x(0) on the boundary 
of X . This results in evolving the boundary of X . More precisely, let S 
be the boundary of Xq, and X(t) = ^oUo<Kt{^( x 0'*) • x o "= ^o}- Then 
X(t) = Reach q (X , [0, r]) (see Appendix A for a proof). This is the idea un- 
derlying the computational approach, called face lifting method, described in Q . 
In this section, we shall study this in considerably more general setting, and de- 
scribe an algorithm, called generalized face lifting method. It may be noted that 
the method described in || does not assume that a global invarience require- 
ment is imposed; hence the algorithm of || computes only Reach q (XQ, [0, r]). 
However, we shall extend our method for computing Reach q (XQ, [0, r]) for com- 
puting an approximation to the set Reach q (Xo), where the approximation can 
chosen to be either under-approximation or over-approximation. 

Further, let Sq C So be the set of those boundary points of Xq, such that 
the solution with initial point in Sq extends into Xq = R™\Ao , i.e., 

S£ = {x G S Q : 3e = e(x ) > such that g (x o , t) G X%, Vt G (0, e)}, (7) 

and define X + (t) = Ao lJ 0<t<T {(/' g (xo, t) : xo G Sq}. It turns out that 
X + (t) = X(r) = Reach q (X , [0, r]) (see Appendix A). Therefore, in order to 
find Reach q (XQ, [0, t]), we have to find only those solutions for which the initial 
conditions are in Sq . 

If Sq can be found explicitly (by inspection of f and Sq), then the problem 
reduces justifiably to finding the solutions with initial conditions in Sq . Other- 
wise, we may have to find a means of obtaining an outer approximation to Sq , 
i.e., a set Sf satisfying Sq C S^ C Sq. We suggest a way to find such a set. For 
this purpose, we assume that Xq is specified as follows: there is a continuously 

o 

differentiable function, I : R" — * R, such that if x G Xq then £{x) < (where 

o 

Ao denotes the interior of Ao, i.e., the largest open set contained in Xq), and if 
x G Xq then £(x) > 0. Hence £(x) =0 on So- Define 

S+ = {xgS : V^(x)./(9,x)>0}. 
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It may be easily shown that Sq C (see Appendix A). Futher, we observe 
that the set of points for which W(x) • f(q,x) > 0, constitutes an inner ap- 
proximation to Sq. The method of finding an as described above can be 
extended easily to the situation where Xq is specified as the set of intersection 

o 

of finitely many sets of the form £j(x) < 0, with the strict inequality in X . 

With this notation, we proceed to describe a schematic algorithm for finding 
X + (t) for t > 0. In order to exploit the semigroup property of the reach set, 
i.e., Reach q (Xo, [0,t]) — Reach q (Reach q (Xo, [0, s]), [0,t — s]), for any t and s, 
with < s < t, we consider a sequence of time intervals [0, n], [ti,T2], [t2,t 3 ], 
etc., where < t± < T2 < T3 < . . ., and Tj — > 00, as i ^ 00. It may be con- 
venient to choose for some r > 0, Tj = it, i = 0, 1,2,.. ., although we do not 
require such an assumption in the algorithm. 



Procedure for Reach q (X , [0, r]) 

initialize: X+(0) := X , F Q := and i := 
while t > n and Fj ^ do 
if r < Tj+i then Aj := r — r» else Aj := r,+i — Tj end if 
:={0 a (x,t) : 46 [0,AJ, x e FJ 

/* this computational step requires special attention! */ 
X+(Ti + Ai):=X+(Ti)\JTi 
F i+1 := {0(x, Aj) : x e F}\A+(r 4 ) 
i := i + 1 
end while 



In the above schematic, the initialization step "Fo := " could be replaced 
with "F := as it is not necessary to initialize F to 5+ . However, before 
proceeding further, it must be mentioned that with reference to this algorithm, 
we assume that the computational steps "7j := {</>q(x, t) : t E [0, Aj], x e F}" 
and "Fj + i := {0(x, Aj) : x e Fj}\A+(rj)" can be performed effectively. 

We now describe a method for extending this procedure to another pro- 
cedure that computes an approximation to Reach' q {Xf), X qi [0, 00)). But since 
Reach q (Xv) = Reach' (X ,G(q), [0, 00)), the procedure to be described below 
finds an approximation to Reach q {X ), when X q = G(q). The approximation 
can be chosen to be either under-approximation or over-approximation, deped- 
ing on a flag "under .approximate" , passed as input to the algorithm (refer to 
the schematic algorithm described below). 

To facilitate the discussion, we consider the flow tp q (x, t) , t > 0, x e X = R", 
corresponding to the differential equation 

f = -/(,. x) (8) 

with the initial conditions x(0) = xq € X = R". Now since the flow ip q 
has the opposite direction to that of <j> q , for any t > and two subsets X\ 
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and X 2 of E", we have, X 2 = 4> q (X u t) if and only if X x = tp q (X 2 ,t). The 
function ij) q will be used in the procedure for computing an approximation to 
Reach' '(X , X q , [0, oo)), for finding, in each iteration, the set of initial conditions 
in Fi corresponding to which the solutions for the time interval [0, Ai] may vi- 
olate the global invariance requirement. With this, we present the schematic 
algorithm as follows: 



Procedure for Reach' q {Xq, X q , [0, oo)) 

precondition: Xo C X q 

boolean input flag: under -approximate 

/* under ^approximate = 1, if under-approximation is preferred */ 
/* by default, the procedure over-approximates the reach set */ 

initialize: X+(0) := X , F Q := S+ and i := 
while Fi ^ do 

A 4 := Tj+i - n 

Ti := {<f> q (x,t) : t€ [0,Ai], xef,} 

if Ti C X g then /* global invariance not violated */ 

X+( Ti + Ai):=X+(Ti)\JTi 

F t+l :={^(x,A0 : x G Fj\X+(r,) 
else /* global invariance violated by at least one trajectory */ 

Ui := T\X q 

V t :={^,(x,t) : (6 [O.AJ, x e C/J 
7V:={0,(x,t) : te [O.AJ, x e F/} 

if under .approximate = 1 then 

X+fo + AO :=X+(rOU^ 
gIsg 

X+(r 4 + A t ) ^X+^OUTi 
end if 

:={^(x,A0 : x e F/}\X+( Ti ) 

end if 

i := i + 1 
end while 



As with the previous algorithm, the initialization step "Fo := S^" could 
be replaced with "F := Sq" , depending on convenience. Termination of this 
procedure may be guaranteed, if certain assumptions are satisfied. The required 
assumptions are as follows: 

1. the set Y q = X q is compact; 
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2. X is closed in K n ; and 

3. for every Xo £ Xq, Y q does not contain any w-limit points of the trajectory 
of (j> q (xQ,t), t > 0, more precisely, 

9 (x o ,oo)P|r g = 0, xo e i . 



Under these assumptions, it can be shown that (see Appendix C) there is a 
T m ax > 0, such that for every x £ Xo, there is a r(x) > with < t(x) < r TOQ2; 
and 0g(x, t(x)) X g . Obviously, with such a r max , we have r x < r max , for 
each x £ Xq, where r x is as in (0). Termination of the algorithm may be de- 

(i) 

duced from this as follows: Consider a sequence of sets F , i > 0, defined 
as Fq 1 -* = {x £ S*o : </> 9 (x, Tj) € Fj}. It is easy to check that i*g is a non- 
increasing sequence of sets, i.e., Fq C Fq . Let x £ Xq, and let k be a 
nonncgative integer such that Tk < r x < t>+i. Now, during the /c'th iteration 
of the while-loop, the global invariance condition (the condition that Tk C X q ) 
is violated, since <p q (x,t) £" X q , for some t such that r x < t < r^+i. Let 
r x ^ [ T x, T k+i) be a time instant such that 0fc(x, r x ) ^ Xq, and let i x = — t^. 
With this choice of i x , we have both < t x < and </> g (x, Tfc + t x ) g" 
Thus, g (x, ^ + t x ) £ Uk = Tk\X q (refer to the schematic above). Let 
y = (j) q (x, r fe + i x ), so ? (x, r fe ) = ip q (y,t x ). But, since y G C4, ^ q {y,t x ) £ V k 
(refer to the schematic). Hence </> g (x, Tfc) e Vfc and </> g (x, r^) ^ F' k = Fk\Vk- 
Therefore </> g (x, Tfc + i) Ffc + i, implying that x g" F By the monotonicity 

of the sequence of sets Fq 1 \ we have x ^ Fq 1 \ for any i > k + 1. Finally, let 
m > be an integer such that r m < T max < r m+ i. Now since r x < r max , 
x ^ F(j m+1 \ for any x £ X . Hence F^ n+1 ^ — 0, so F TO+ i = 0, satisfying 
the terminating condition of the while-loop of the algorithm, after at most m 
number of iterations. 



In the rest of the paper, we discuss approximation of the set Ti appearing in 
the second line of the while-loop of the above schematic (see also the schematic 
for Reach(X , [0,r])), when f(q, x) is of the form /(g,x) = A q x, where A q is a 
constant n x n matrix with real entries and the initial set Xq is a polyhedron. 
For convenience, we drop the subscript q in A q , and deal with the linear system 

*c(t) 
dt 

where A is a constant n x n matrix with real valued entries. 



4 Representation of the Reach Sets 
4.1 Preliminary Discussion 

In this section, we discuss representation of the reach set Reach q (Xo, [0, r]), 
such that boolean set operations such as union and complementation can be 
performed efficiently. The most convenient representation schemes appear to be 
representation in terms of the following classes of subsets of M. n : 



14 



1. Polyhedral Sets : These are the sets which can be written as union of 
finitely many polyhedra. A polyhedron in K." is a set which may be 
expressed as the intersection of finitely many closed half spaces Jl8|], i.e., 
a finite intersection of sets of the form {x € R™ : a T x < b}, where a £ R™ 
and fi e 1. It may be observed this class of sets corresponds to (and 
includes) the class of definable sets in the theory of linear inequalities of 
R, i.e., the theory obtained when and 1 are the constant elements, + 
and — are the (binary) function symbols, and < is the (binary) relation 
symbol (see Appendix A). 

2. Subalgebraic Sets : These are the sets which can be written as union of 
finite number of sets defined by polynomial inequalities, i.e., sets that can 
be expressed as union of finitely many sets, each of which is the intersection 
of finitely many sets the form {x S R™ : p(x) < c}, where p : R n — > R is a 
polynomial function with real-valued coefficients and c € R. This class of 
sets corresponds to (and includes) the class of sets that are definable in the 
theory of R when viewed as an ordered field. The theory of ordered field 
of R is the theory obtained by extending the theory of linear inequalities 
of R by including a binary function symbol for representing the product 
of two real numbers, denoted by • (see Appendix A). 

It may be observed that, for a general flow function ip q , the reach set 
Reach q (X , [0, r]) may not be exactly representable in any one of these classes, 
even if Xq is. Hence, as an alternative, we may have to settle for an approx- 
imation of the reach set by sets that belong to the respective classes. In this 
context, we shall restrict our attention to over-approximation of the reach set 
with polyhedral sets. 

Specifically, we describe an algorithm for reach set over-approximation with 
polyhedral sets, when the system dynamics (within a discrete state q) are spec- 
ified as 

— -Ax, x(0)el , (9) 

where the initial set Xq is a polyhedron and A is a constant real valued n x n 
matrix. Based on the discussion of the last section, we consider the problem of 
approximating the reach set of the solutions of the equation 

— = Ax, x(0)GF , (10) 

where Fq is a face of Xq. 

The main references on approximate reach set computation using polyhedral 
sets appear to be Q and |lj] . In |l(J , the initial set is assumed to be convex (not 
necessarily a polyhedron), and for each time instance, each slice of the "reach 
tube" (i.e., for each t g [n, Tj+i], the set (j){Fi,t)) is approximated by polyhedra. 
The method described in |l^] assumes Fi to be a polyhedron, and approximates 
the tube for the time interval [r i; r i+1 ], i.e., the set Tj = U 0<t<r +1 _ r . 4>{Fi,t), by 
first constructing a polyhedral over-approximation of the reach tube for the time 
interval and further over-approximating the resulting polyhedron by a "griddy" 
polyhedron, i.e., a set that may be expressed as a union of unit hypercubes with 
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integer left-most vertices 
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Being a convex set, any polyhedron that over-approximates the reach set 
contains the convex hull of Tj. An algorithm based on this observation (as 
given in |l^]) is to find a "bloated" convex hull of Fi and Fi+\. This method 
seems to be simple and straightforward, and gives reasonable approximations in 
a rather short time (refer to fl^l ). But we may note that those faces of the con- 
vex hull of Fi and Fi+i, for which the solution set Tj lies on one side, need not 
be bloated. By shrinking the convex hull, an under-approximation is obtained 
by the same method. 

If the time steps are constant, i.e., = kr, for some r > 0, as in then we 
have to find approximation Pq only for Tq, since symbolically, T k = e^ AT 'T)-—i, 
for k > 1, and therefore, T k = e^ kA ^T = e^ AT ^T . Thus, P k = e^^Po is the 
required approximation, and if X 3 ■ , j = 1, 2, . . . , TV, are the normals of the faces 
Pq, then the normals of the faces of P k are given by e^ ATTk ^ Xj. From this, we 
may infer that in the subsequent iterations of the algorithm of ]l2] | , the convex 
hull of Fi and Fi + i need not be computed. 



4.2 Over-approximating the Reach Set with Polyhedra 

In this section, we consider the problem of over-approximating the reach sets, 
when the dynamics are specified as 

^ = Mt), t > o, 

with the initial conditions x(0) £ Xq. The solution is given explicitly by x(i) = 
e( At ^x(0). We assume that the initial set Xq is a polyhedron, defined by 

\Jx-hi<0, i=l,2,...,K, (11) 

where k is a fixed positive integer. We assume that Xi £ M™, |Ai| = 1 and hi £ R, 
for 1 < i < k. Then the boundary set of X Q consists of sets Cj, 1 < j < k, of 
the form 

Ajx — hj — 0, and 

A^x — hi < 0, for 1 < i < n and i =/= j . 

Further, we consider only those polyhedral initial sets Xq and matrices A that 
satisfy the following assumptions: 

{A\) X is compact. 

(^2) The set Sq , as defined in ([?]) corresponding to the flow function <p q (x, t) — 
e^^x, can be expressed as Sq — Cj 1 \J Cj 2 U ■ ■ • U Cj m > where 1 < j% < 
32 < ■ ■ ■ < jm < such that for each i with 1 < i < m, Xj.Ax — hj t > 
Si > 0, for every x £ Cj t . 

Referring to these assumptions, it may be mentioned that the main limitation 
of this method appears to be the restriction imposed by assumption Ai on the 
initial set Xq and the matrix A. In contrast, the method described in |12 docs 
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not require such an assumption. But for our purposes, we need this assumption. 

Now, let j be an index such that Cj C Sq, as in assumption (.A2). For 
simplicity and definiteness, let us assume that j = n and put Fq = C K . Hence, 
Fq is described by 

Xfx-hi < 0, l<i<K-l, and 
X^x-K = 0. 

After some manipulations and rearrangements, if necessary, the above system 
of linear constraints can be transformed into an equivalent system of linear 
constraints, which also describes Fq but is of the following form: 

afx-bi < 0, l<i<fc-l, and 
a^x - b k = 0, 

where k < n, &k = X K , and for 1 < i < k, |aj| = 1 and afa^ = 0. It may be 
observed that, since |Afe| = 1 initially, we have |a&| = 1 as well. Also, by as- 
sumption (A2), there is a 6 > 0, such that for all the vectors x £ Fq, Ax > 5. 
Further, we assume that the system of linear constraints in ( |l2| ) is consistent 
and that each of the constraints is linearly independent from the remaining con- 
straints. 

With this construction, our objective is to describe an algorithm to find a 
polyhedron Pq as an over-approximation to the set 

T = {e^x : x £ F , ie[0,A]}, (13) 

where A is a small positive number for which the the following condition holds: 

(Ci) for some S > 0, Ae^ At "> x > So, for all x S F and for all t £ [-A, A]. 

The existence of such a A can be assured by our assumptions {Ai) and (A2)- 
Later, we shall derive an estimate for such a A, for any fixed So > such 
that So < S. In comparison, the method of Q does not impose any such 
restrictions on A. But, it may be mentioned that, in both cases, the accuracy 
of approximation of either method depends on how small a value is chosen for 
the parameter A. The two conditions (C2) and (C3) stated below follow from 
condition (Ci): 

(C 2 ) a^e^'^xo - b k > 0, for all x £ F and for all t such that < t < A; and 

(C 3 ) a^e^xo - b k < 0, for all x £ F and for all t such that -A < t < 0. 

Before proceeding to describe our algorithm, we note that the set -Fa = 
{e^ A ^xo : xo £ Fq} may be described by the system of linear constraints 

af(A)x-6. t < 0, i = l,2,...,k-l, 
a£(A)x-fe fc = 0. 

where a^(A) = e^~ A A ^a.i, 1 < i < k. Dividing throughout the last equation in 
the above system by |afe(A)|, we get 

af(A)x-6 l < 0, i= 1,2,...,*- 1, 
b!x-6' fe - 0, 
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where bt = , afc ^|| and b', — , b ,\-,, . Now, we observe that, since aT A ) i s 

K K(A)| K K(A)| ' ' 

invertible and |a,| = 1, each a^(A) is nonzero, and since a.k is orthogonal to 
aj, for 1 < i < k — 1, afc(A) and a.;(A) are pair-wise independent. So, after 
subtracting from each of the inequalities of the above system an appropriate 
constant times the last equation, and after normalization, the above system 
of linear constraints can be transformed into the following system of linear 
constraints: 

hfx-br < 0, i = l,2,...,fc-l, 

where, for 1 < i < k, |bj| = 1, and for 1 < i < k—1, hfhk — 0. Now, by the con- 
dition (C 3 ), for < t < A and x G F , a£(A)e( A *) -b k = a?V A (*- A » -b k < 0. 
Hence, we also have, for < t < A and x £ F , bje^^xo — b' k < 0. 

In this notation, we now describe a schematic algorithm, the output of which 
is a set of 4k many parameters, consisting of pairs of vectors and constants, 
(/ii, Ci) and (yi, di), i = 1,2,..., 2k, where for each i, Vj, € R n and Cj, d{ 6 R, 
such that the polyhedron Po of intersection of the 4k half-spaces defined by 
Li = {x S R" : /zf x - c, < 0} and L- = {x e R™ : i/f x - ck < 0} is an 
over-approximating polyhedron for the set T as defined in (|l^) . 




18 



Schematic Algorithm for Finding Over-approximating Polyhedron 

Reach Set : T a = {e (At *>x : t e [0, A], x € F } (refer to @) 
Output : 4fc vectors /tj, i/j and 4fc real numbers Cj, dj, 1 < i < 2k 

/* the 4fc half-spaces are Lj and L' lt I < i < 2k, where */ 

/* Li = {x € E n : ??i(x) = /if x - a < 0}, and */ 
/* L^jxe E" : ^(x) = i/fx - di < 0} */ 

for i = 1, 2, . . . , k — 1, do /* to find (/tj, Cj) and (z/j, dj) */ 
U := inf{l : (af x - 6,) - Z(af x - b k ) < 0, Vx e T }; 

Mi ; — a i — ^i a fe! Cj - = bi — libk] 

l[ := inf{f : (bf x - 60 + Z'(bf x - b' k ) < 0, Vx G T }; 
i/i := b, + l'^; d, = b\ + l^b' k ; 
end for 

fi k := -af ; c fc = -b k ; 



:= bf ; d k = b' k : 



l k := inf{Z : (af x - b k ) -K0,VxeT o }; 

Mfc := a fc^ c k+l ;= ^A; + Zfe! 

l' k := inf{/ : -(bf x - b' k ) - I 0, Vx e T }; 
^fe := — bfc; := — b fe + i fe ; 

for i = 1, 2, . . . , k — 1, do 

/* to find (fx i+ k+i,c i+k+ i) and (i/i+k+i, d i+ k+i) */ 
l i+k = inf{Z : (af x - b t ) - I < 0, Vx e T }; 
Mj+fc = a^; = bi + li+ k \ 

l' i+k = m£{V : (bf x - b'i) - V < 0, Vx e F }; 
Vi+k = b^; rfi+fe+i = b' L + l' i+k ; 
end for 



We first observe that both U > — oo and l\ > — oo. If af Ax > 5 > 0, as 



in assumption A2, then both It < 00 and l\ < 00, as we shall see in Sec. 4.3 
below. If A is small, such that the set Tq is nearly a polyhedron, then the 
above method gives reasonable results. Referring to Step 1 of the first for-loop 
of the above algorithm, the reason for choosing Z, to be the infimum over all 
I for which (af x — bi) — Z(af x — b k ) < 0, Vx e T , is obvious: if Pi and P2 
are two polyhedra such that Pi is obtained by the above algorithm and Pi is 
obtained by replacing a constraint /if x — a < 0, where /t, and Ci are as in Step 
2, with another constraint Af x— ft,; < 0, where A.; = a,; — Za^ and hi = bi — Ibk, 
for some I > h, then Pi C P2-0 (An analogous statement holds for each of the 
remaining parameters, l[, l k , l' k , l k +i and l k+i , as in the algorithm.) 

3 This follows from the observation that if both (af x — bi) — /(af x — ftj.) < and — (af x — 
bk) = ^f x ~ c k ^ 0) where I S R and x S R n , x not necessarily in To, then, for any h > I, 
(af x - b t ) - h(a£x - 6 fe ) < 0. 
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It can be shown that the polyhedron included in the k + l half-spaces, speci- 
fied by Li, L27 ■ ■ ■ j Lk, Lk+i, is a bounded polyhedron (see Appendix D). We may 
also note that the hyperplanes /ifx — Ci = and vjx — di = 0, 1 < i < k— 1, are 
obtained by rotating the hyperplanes a^x — bi = and b^x — b\ = about their 
corresponding intersection with the hyperplanes a^x — bk = and b^x b' k = 0, 
respectively; whereas the hyperplanes /i^ +1 x — Ck+i = and v^ +1 x — dk+i = 
are obtained by translating the hyperplanes a^x — bk = and b^x — bk = 0, 
respectively. Similarly, for i = 1, 2, . . . , k — 1, the hyperplanes ^ +i x ~ c k+i = 
and — dk+i = are translations of the hyperplanes a^x — bi = and 

hf x — = 0, respectively. 



In remaining part of the section (in Sec. O below), we shall derive upper 
bounds for the numbers U and l[, 1 < i < 2fc — 1, that are defined in the 
schematic. If all the U and are replaced with their corresponding upper 
bounds, lj and l[ respectively, then we obtain conservative estimates — Aii Cj, £j 
and respectively - for /Ltj, Cj, i/j and d^. After this replacement, we obtain a 
modified algorithm for over-approximation with polyhedra. It may be observed 
that these estimates overcome the difficulties in finding the infima required 
for obtaining Zj's and I'^s. But unfortuna tely , the conservative upper bounds 
that we derive for Zj's and Z-'s (as in Sec. |4~3| ) may turn out to be very large. 
However, better accuracy may be obtained if the modified alogrithm is used 
in conjunction with that of jl^] (see Fig. || in Sec. 5). More precisely, the 
intersection of the polyhedron obtained by the method described here with that 
obtained by the method of p2[ gives a smaller over-approximating polyhedron, 
as will be discussed in Sec. 5, while illustrating these algorithms with simple 
examples. 

4.3 Upper Bounds for U and /' 

In this section, we derive upper bounds for the numbers Zj and l[ appearing in 
the schematic algorithm. We begin with the following result: 

Claim 1. For i = 1, 2, . . . , k, with a^ and bi as in ( |l2] ) and xo G Fq and 

t £ (0, A], we have 



(afe^xo - a[ x ) 



t 

where M = max xoeFo {||x ||}. 



< Mo||a^4||e<IWI*>, 



To prove the claim, we first note that (af e^ At ^XQ — afxo) = t x af Ae( Ae ^x.o, 
for some 9 with < 9 < t, where may depend on t. Hence, we have 



(aJee iAt) Xo - afx 



= \af Ae {Ae) x \, where 9 is such that < 9 < t 

< \afAe^\ 

< M \\afA\\e 



(MIIA) 



With hi and b[, we have an analogous result, but with a slight modification, as 
in the following: 
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Claim 2. For i = 1,2, ... ,k, with and b\ as in ( p^ ) and xo G Fq and 

t € [0, A), we have 



(bf e^xo ~ bf e^ A )x ) 



(A-t) 

where M = max XoeFo {||x ||}. 



<M \\b?A\\e 



(MIIA) 



We have (bfe^Xo - bf e( AA >x ) = (t - A) X bf Ae^x , for some 6> with 
t < 6 < A. Hence 



(bf e(^)xo - bf e( AA )x ) 



(A-t) 



= |b^e (A9) x |, for some with t < 8 < A 

< ||bf Ae^H X ||xo|| 

< M ||bfA|| e (ll A ll A ). 



Finally, we shall make the following claim, before deriving upper bounds for 
h's and Z-'s: 

Claim 3. Assume, for some S > 0, af Ax > S, for every x G F , 
and let So be such that < 5q < S. Then there is a A > such that 
af Ae^xo > S , for all x G F and for all t G [-A, A]. 

In order to prove the claim, we shall derive a conservative estimate for a A > 
for which the claim holds. Fix a So such that < So < S. We have 



TAe^ 



x 



af Ax 



af A 2 e (As) x ds, where x G Fq, 



r\t\ 



> afAxo- / ||af||P||V^ll s )||x ||d S 

Jo 

> afAxo-Moll^KeWI'l-l), 
where Mo = max Xo gi? {||xo||}. Therefore, if A is such that 

M \\A\\(e^ A - 1) < S ~ So, 

then the claim holds. 

Upper Bounds for U and l[, 1 < i < k — 1. Recall that, for 1 < i < k — 1, 

h = mf{l : (af x - h) - /(af x - b k ) < 0, Vx G T }. (15) 

Referring to the right hand side of (|l~5|), for every x G F , since af x — bk = 0, 
for any real number /, we have (af x — bi) — /(af x — bk) — (af x — hi) < 0. 
Therefore, we may assume x G T \F , and (af x — 6 fe ) > 0. Thus li has to be 
chosen such that 

sup 



u 



Rewriting the above, we have 



<£T \f (afx-6 fe ) 



l = gu (af e^xo - 6Q 
o<t<A (af e( A ')x - b k ) ' 
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Dividing the numerator and denominator by t, this may be written as 

ai(x ,t) 
li = SU P 1T< IT' 

0<t<A PH x O) t) 
xoG-Fo 

where ai(x ,t) and (3i(x ,t) are given by 
(af eMxo - h) 



Q!i(xo,i) = 



t 



and A(x ,t) = /?(x ,t) = t^ 6 ^* ^ . 



But, since af x — h <0, for every x € F , we have 



(afe^)x -afx ) 



(afe^xp-b,) < (ajV^Xp-afxp) < 
t ~ t 



By Claim 1, we have 

sup ai(x ,t) < M \\afA\\e^ A \ 

0<t<A 

As to the denominator, since a^x — — 0, for every x 6 F 0l we have 



(&W At) xo - M = (^e^xo - a^x ) = ^^{Ae)^ > § 
t 



T 



where, 9 is some number in the interval (0,i), and the last inequality is due to 
Claim 3. Hence 

inf /3(x ,t) > (5 > 0. 

0<*<A 

Combining both inequalities, we have 

h < su Po<t < A: Xo£Foaa (x ,t) < . _ MoHafAHedl- 4 ^) 
info<t<A, x e_F /?i( x o, i) 

Likewise, /' must be chosen such that 



, a t (x ,t 

k = S UP -577 7T • 

o<t<A Pi(, x o,rJ 
xoefo 



where a^(x ,t) and /3-(x ,t) are given by 



o / < (x ,t)= (b ^? X ° ^ and #(xo,t) =/?(*>,*)- l " ; 



bZV^xo - &' fe ) 



(A-t) 

Now, since bf e ( AA )x - 6- < 0, for x e F , 

(bf e( A *)xo - ^) < (bf e( A *)x - bf e( AA )x ) 



(A-t) 



(A-t) 



< 



(A-t) 
(bfe( A *)xo-bfe( AA )) 



(A-t) 

< M ||bf A||e(ll A ll A ) , by Claim 2. 
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For the denominator function, since b^e^ A 'xo — b' k = 0, for every x G F , we 
have 

(b^ (At) x - b' k ) = (b^)xp - b^ A >x ) _ _ h T Aei A0 )x 
(A-t) {A-t) k °' 

where is some number in the interval (t, A). Now, 

a T„(-AA) AjAd) v s 



.aje(-^)|| - ||a^e(-^)||- 
Hence, we have the following upper bound for l[: 

t < su Po<t < A , XoeFo a' t ( Xo ,t) < f , = MollafAHe^ll^ 

~ info<t<A, xoG-Fo $( x 0, t) ~ 1 Si 

Thus, we obtain the following upper bounds: 

l\ < l\ = M ° l|a * 7 e , l<t<*r-l, where 

01 

t = s ° 
1 lla^-^H' 

Upper Bounds for Ik, l' k , h+i and l' k+i , 1 < i < k—1. We have k + k — inf{l : 
(af x — 6i) — Z < 0, Vx e T }, and Z fe = inf{Z : (a^x - 6 fc ) - I < 0, Vx e T }. 
Let ?2fe = ^fe, so we may consider Zfc + j, for 1 < i < k. We have to choose lk+i 
such that 

Z fe+4 = sup (af x - bi) = sup (af e (At) x - 

xeT 0<t<A 
x eF 

Now, since (af x — bi) < 0, for every x € F , we have 

sup (af e (At) x - 6j) < sup (af e (At) x - af x ) 

0<t<A 0<*<A 
xoSfb x eF 

(af e( A *)x - af x ) 



< sup t x 



0<t<A t 
x eF 

(af e (^)xo - af x ) 



< Ax sup 



0<t<A t 

< AxM ||a^||e(ll A ll A ). 
Similarly, for 1 < i < k — 1, we have to choose as follows 
4 +l = sup (hf x - 60 = sup (bf e^x - 

xeT 0<t<A 

x eF 
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A calculation similar to the above shows 

sup (bf e^Ox - 6<) < AxM ||bfA|| e (ll A ll A ). 

0<t<A 

Finally, for l' k , we have 

4 = sup [-(b£x - b' k )] = sup (b£e^ A >xo - b£e^>xo) 

xST 0<t<A 
xoG-Fo 

Another sequence of similar calculations shows 

sup (b^ A >xo-b£e^>x ) < AxMoHb^lledl- 4 ^). 

0<t<A 

So, to collect all the estimates, we have 

l k < l k = M a A\\ a T k A\\e^ A ^ 

l' k < i' k = M A\\h T k A\\eJ^ 

lk+i < lk+i = Af A||af A\\eS^ A \ l<i<k-l 

l' k+l < l' k+i = M„A||bf A\\e^ A \ l<i<k-l 



5 Illustration 

5.1 Example 1 

We first illustrate the schematic algorithm presented in Sec. 3 with the help of 
an example taken from |^|. Let 

X . 2 a,b>0, (x(0),y(Q))6X = [0,l]x[0,l]. 

So, for t > 0, x(t) = x(0) + at and y(t) = y(0) + bt. It is easy to check that 
So = {!} x [0> !] U [0, 1] x {!}■ Therefore 

X+(t)=X (J {{l+at,c+bt) : < c < 1} (J {{c+at, 1+bt) : < c < 1}. 

0<t<r 0<i<r 

This is illustrated in Fig. 1. 

5.2 Example 2 

In this example, we illustrate the over-approximation algorithm. Let x(t) = 
[x 1 (t),x 2 (t)] T £ K 2 satisfy 

f I ~ X2 \ (x 1 (0),x 2 (0))GF = [l,V2]x{0}. 

The solution is given explicitly by 







cost 


— sint 




' xi(0) 


. x 2 {t) _ 




sini 


cost 




. x 2 (0) 



24 




Figure 1: Illustration of the Generalized Face Lifting Method: X + (t) 
X {jA{jB. 



With A = 7r/6, for the time interval [0, 7r/6], the solution set in a parametric 
form is To — {(acost, asint) : a £ [1, \/2], < t < 7r/6}. For this example, in 
the notation of Sec. 4, we have 





' 1 " 




' -1 " 




' " 


ai = 





, a 2 = 





, a 3 = 


1 



y/2, b 2 



-1, and b?, 



0, 



and 
bi = 



r Vs - 




r n/3 - 




' -0.5 " 


2 

0.5 _ 


, b 2 = 


2 

-0.5 _ 


, b 3 = 


L 2 J 



, b[ = a/2, b' 2 = -1, and b' 3 = 0. 



The set To and the result of the algorithm with li and l[ exactly found as in the 
schematic of Sec. 4 are shown in Fig. 2, where the horizontal axis corresponds 
to x\ and the vertical axis corresponds to x^. 



A 



A 




-> 




-> 



Figure 2: Illustration of over-approximation algorithm of Sec. 4 with Example 
2. Left: The set To. Right: Result obtained by the schematic algorithm of 
Sec. 4, with exact values for U and l[. It may be observed that only 5 of the 
4fc = 4x3 = 12 hyperplanes that the the algorithm computes are required for 
this example, provided that and U can be found exactly. 



Now, recalling the notation of Sec. 4, we may take S = 1, Sq 



2 ' 



It is 
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easy to check that M = y/2 and ||af A\\ = \\hjA\\ = 1, for 1 < i < 3. Hence 
y/2 x 1 x e^/ 6 2a/2 x e 7 ^ 6 



\/3/2 



V3 



2.7566424 



and the same upper bound holds for since <5i = <5o in this example. For 
and I'u we have the following bounds: 

h = l'k < h = V2 x | x e^ 6 w 1.249999 

In this example, the same upper bound, given by Zfc, holds for and 
Combining all this, the polyhedron is the intersection of the half-spaces given 
by 7?i(x) < and f?,-(x) < 0, 1 < i < 6, where 

- 0.5122958 x x + 2.8873223 x 2 

- 2.2443466 x x + 1.8873223 x 2 

- Xi + VSx 2 , 

/3 

0.5xi ~ — x 2 - 1.249999. 



?7i(x) = x x - 2.7566424 x 2 - y/2, 


^(x) = 


77 2 (x) = - xi - 2.7566424 x 2 + 1.0, 


= 


%(x) = -x 2 , 


^s(x) = 


j? 4 (x) = x 2 - 1.249999, 


^(x) = 


775 (x) = xi-V2- 1.249999, 




?? 6 (x) = - xi + 0.249999, 


^e(x) = 



73 



;i +0.5x2 - V2- 1.249999, 



-^-x x - 0.5x 2 + 0.249999. 



In this example, ^(x) < 0, f76(x) < and r/g(x) < 0, are redundant. 

Fig. H illustrates the result of the algorithm with the first 8 half-spaces de- 
fined by ??i(x) < and ^(x) < 0, 1 < i < 4, where the dashed lines correspond to 
the lines rji(x) = and 77-(x) = 0, and the polyhedron included in their intersec- 
tion is shown in thick lines. The point of intersection of r/i(x) = and 774 (x) = 
is (4.8600138, 1.249999), of r/ 4 (x) = and r?i(x) = is (4.2845099, 1.249999), 
and of 77 2 (x) = and r? 2 (x) = is (0.575162, 0.154114). So the vertices 
in the counter-clockwise order are given by (V2, 0), (4.8600138, 1.249999), 
(4.2845099,1.249999), (V3/V2, 1/V2), (V3/2, 0.5), (0.575162, 0.154114) and 
(1,0). 

Fig. U shows the result when the remaining two half-spaces defined by 
%(x) < and 775 (x) < are also used. Finally, Fig. || shows the intersection of 
the polyhedron obtained by the method reported in this paper with that which 
may possibly be obtained by the method of [^2) , the latter being shown in dotted 
lines. The polyhedron corresponding to the method of |l2| is obtained as follows: 
the convex hull of the sets [1 , V2] x {0} and { (a cos ^, a sin |) : o 6 [1, v2] } is the 
polygon with vertices in the counter-clockwise order (-\/2,0), (V3/V%, 1/V2), 
(a/3/2,1/2) and (1,0). Therefore the poly gon is the interesction of the half- 
spaces Ci(x) < and Ci( x ) < 0, i = 1,2, where 

Ci(x)=-x 2 , Ci(x) = 

1, ( V3 \ 

x 2 . 



C 2 (x) = ±( Xl - V2) - ^ - V2j x 2 , C 2 (x 
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Now if e is an upper bound for the bloating parameter e, then the half-spaces 
of the over-approximating polyhedron corresponding to the method of [l2| are 
given by Ci( x ) ^ an( i Ci( x ) < 0, z = 1, 2, where 

Ci(x) = -x 2 + e, 

Ci(x) = -x 1 + V3x 2 + e, 

Ca(x) = 0.70710678x1 + 0.18946869x2- 

C 2 (x) = -0.70710678xi-0.1339746x2 



1 + e, and 
- 0.70710678 + e 



The upper bound for the bloating parameter e as given in |0| works out to be 

e < A/ ^II A II A )-1-||A||A-^P|| 2 A' 

3^ 2 



72: 



,(t/6) 



7T 

1 - - 



6 8 x 36, 
= \/2 x 0.06168464 « 0.087235255 < 0.09. 

Fig. ^ shows the results when an upper bound for the bloating parameter is 
chosen to be e w 0.2, where the dotted lines correspond to the line Ci(x) = 0, 
Cl(x) = 0, C2( x ) = and C 2 (x) = 0. As may be expected, the polyhedron 
of intersection of the two polyhedra - one polyhedron bounded by the dashed 
lines corresponding to the method described here and another bounded by the 
dotted lines corresponding to the method of jl^] - gives better results of over- 
approximation. 




Figure 3: Illustration of polyhedral over-approximation for the soultion set To 
of Example 2: without additional hyperplanes. 



6 Discussion and Conclusion 

An important issue of hybrid systems appears to be computability of the reach 
sets of the continuous variables. From a computational point, both the problems 
of computation and efficient representation of the reach sets of the continuous 
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A 




Figure 4: Illustration of polyhedral over-approximation for for the reach set Tq 
of Example 2: with additional hyperplanes. 

A 




Figure 5: Illustration of polyhedral over-approximation for the solution set 
To of Example 2: intersection of the polyhedron that may be obtained by the 
method of with that obtained by the method reported here. 



varibles are difficult, in general, owing to the limitations of quantifier elimination 
method. In this context, approximation of the reach sets by more convenient 
sets, such as polyhedra and subalgebraic sets, is discussed in the literature. 
In this paper, along with a method for finding the reach sets, an algorithm 
for over-approximation of the reach sets with polyhedra when dynamics of the 
continuous variables are specified by linear differential equations and the inital 
set is a polyhedron. A practical version of the over-approximation algorithm 
is also discussed in this paper. However, it seems that better results of over- 
approximation may be obtained by taking the intersection of the polyhedron 
obtained by the method reported here with that obtained by the method given 
in fl^| . It is hoped that the over-approximation method presented here may be 
extended to systems with more general dynamics and initial sets. 
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Appendix A. Decidable Classes of Subsets of Q x 

W 1 

We introduce here model theoretic concepts for formalizing the notion of decid- 
able class C of subsets of Q x R™ that we were interested in Sec. ??. Since Q is 
a finite set, for a decidable class C of subsets of Q x R n , we may take, for each 
q £ Q, a decidable class of subsets of R™, S q , q £ Q, and choose C to be the 
product class U 9 gq{{'?J' x S q : S q £ S q }. Therefore, we restrict our attention to 
the discussion of classes of subsets of R™. Most parts of this section are taken 
from H (see also [?] , which in turn cites the references jl^, |l| ) . 

Definition 11 A language is a tuple of three sets, C = (L r ,Lf,L c ), where 

1. L r is a set of relations, 

2. Lf is a set of functions, and 

3. L c are a set of constants. 

Definition 12 A model of a language C = (L r , Lf, L c ) consists of a nonempty 
set S, together with an interpretation of the relations, functions and constants. 

We denote a model by (S, L r , Lf, L c ), where the interpretation is not made 
explicit. In the following, let V = {x, y, z, Xq, x%, . . .} denote a countable set of 
variables. 
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Definition 13 A term of a language, C = (L r ,Lf,L c ), is inductively defined 
as follows: 

1. each variable 6 £ V is a term, 

2. each constant c € L c is a term, and 

3. for an m-ary function g £ Lf, where m > 1, and m terms, 6-1,62, ■ ■ ■ , 6 m , 
9(61, 62, ■ ■ ■ , 6 m ) is a term. 

Definition 14 An atomic formula of a language C is either 

1. 6\ = 62, where 6\ are 62 two terms of C, or 

2. p{6\, . . . , 6 n ), where n > 1 and p <G L r is an n-ary relation. 

Definition 15 A first order formula, or simply a formula, of a language C is 
recursively defined as one of the following: 

1. an atomic formula, or 

2. -i0, where is formula and -1 is logical negation, or 

3. 0i A 02; where 0i and 02 are formulas and A is logical and, or 

4- 3x : or\/x : 0, where is a formula, x is a variable, andB (there exists) 
and V (for all) are quantifiers; in this case, each occurrence of the variable 
x in the formula is called a bound occurrence. 

Definition 16 The occurrence of a variable in a formula is free, if it is not 
bound. A sentence in a language C is a formula with no free variable. A theory 
of C is a subset of sentences. 

For a model S of the language L, we shall be particularly interested in the 
theory defined by the set of all sentences that are true in S. To emphasize this, 
we refer to this theory as the theory of (S, L r ,Lf, L c ). 

Definition 17 Let £ be a language and S be a model of C A set X C S n is 
said to be definable in the language C, if there is an n-ary formula <j>(xi, . . . , x n ) 
such that X can be written as X — {(xi, . . . ,x n ) e S n : <fi(xi, . . . ,x n )}. 

Definition 18 Let C be a language, S be a model of £, and C be the class of 
definable sets. Then C is said to be decidable, if the theory of (S,L r ,Lf,L c ) is 
decidable, i.e., there is a decision procedure that, given an C-sentence 0, decides 
whether belongs to the theory of (S,L r ,Lf, L c ) or not. 

Examples 

1. The theory (R, {<}, {+, — }, {0, 1}) is the theory of linear constraints with 
integer coefficients, denoted by lin(l). The sets defined by these formulas 
are called polyhedral sets. 
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2. The theory (R, {<},{+,—,•}, {0,1}) is the theory of polynomial con- 
straints with integer coefficients, denoted by OF(R). The sets defined 
by these formulas are called subalgebraic sets. 

The result stated below is due to A. Tarski @|: 
Theorem 1 The first order theory OF(R) is decidable. 

Definition 19 Let C be a language and let S be a model of C. We say the 
theory of(S, L r ,Lf, L c ) admits quantifier elimination if every first order formula 
of(S, L r , Lf, L c ) is equivalent to a formula of (S, L r , Lf, L c ) without quantifiers. 

Examples: Decidability and Quantifier Elimination 

1. The theory OF(R) consisting of (R, {<}, {+, -, •}, {0, 1}) admits quanti- 
fier elimination and is also decidable. 

2. Let OF exp (M.) be the theory consisting of (R, {<},{+,—,-, exp}, {0, 1}), 
where exp, representing the exponential function, is a new function sym- 
bol. This theory does not admit quantifier elimination, and it is not known 
whether this theory is decidable. (See [^5|). 

Appendix B: Reach q {X , [0,r]) = X(t) = X + {r) 

We assume that Xq is closed and its boundary, denoted by So, consists of a 
finite union of smooth surfaces. Let 

X(r)=X |J {<Mx ,i) : x G S }. 

0<t<T 

Also let 

S +={xe5 : 3e = e(x) > such that 0(x, t) £ X C , Vi £ (0, e)}, (16) 
and X+(t) = A Uo<t< r {^( x o^) : x o G S$}. 
Proposition 1 X+(r) = X(r) = Reach(X , [0,r]). 

Proof. Note that X + (r) C X(t) C Reach(X o ,[0,T}). So we have to show 
that Reach(X , [0,r]) C X+{r). Let z e Reach(X ,t), t > 0. If z G X Q , 
then z e X + (t). Now suppose z ^ Xq. So there is a Zo € Xo such that 
z = 0(zo,s), for some s with < s < r. If zo € «Sq , then z 6 A + (t). Oth- 
erwise let t' — sup{t : t < s and 4>(z ,t) £ X }, so < r' < s < r. Now 
z' = 0(z o , r') € Therefore, z = <j>{z' , s - t') G A+(r). 

Now suppose that Xo is specified as Xo = {x € R™ : £(x) < 0}, where 

o 

£ : R™ — + R is continuously differentiable. Further, we assume that if x e Ao 

o 

then £(x) < (where Xo denotes the interior of Xo, i.e., the largest open set 
contained in Xo) So, obviously, if x € So then £(x) = 0, and if x £ X§ then 
£(x) > 0. Let 5^ c S Q be defined as 

5+ = {xeS* : V*(x)-/(g,x)>0}. 
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We now show that Sq C S±. To this end, we show that Sq\Si C So\Sq . Let 
x € S Q \Sf, so V£(x) • /(q,x) < 0. Now, since £((f> q (x,0)) = £(x) = 0, and 
since at t = 0, M ^iJ^^ = V£(x) • /(g,x) < 0, we have, for a sufficiently small 

o 

e > 0, £(<f) q (x, t)) < 0, whenever < i < e, which happens only if ^(x, i) 6 Xq, 
Vt e (0, e), implying that x ^ hence x e 5 \5^. Thus S \Sf C 5 \5(|, as 
required. 



Appendix C: Reach' q (Xq, X q , [0, oo)) 

Let X q C K™ with compact closure, and let Y q = X q . Let f : W — > M™ be 
a continuous function defined on an open set VF containing Y q , satisfying a 
Lipschitz condition on W. Let Xo be a closed subset of Y q . We assume that 
for each x G Xo, a function 7 X : R + — > K." exists and satisfies the differential 
equation 

^ = f(7x), t>0, (17) 

with the initial condition 7 X (0) = x. Hence the flow 0(x, t) associated with 
equation (Jl^) is defined for all t > and x <E X . Further, assume that the 
w-limit set of the flow </>(x, t), t > 0, does not intersect F g for any point x E X . 
More precisely, we assume, for x 6 X , 

L a) (x)f|y 9 = 0, (18) 

where 

A„(x) = f|0(x, [i,oo)). 
t>o 

(See [|l6| |l^].) In what follows, we show that if ( |l8| ) holds for every x e X , then 
3r max > (depending on Xq), such that for every x € Xo, there is a t = t(x) 
with < r(x) < r max and </>(x, t(x)) ^ Y" 9 . For this, let < n < t 2 < . . . be an 
increasing sequence such that T k — * oo, as k — > oo (as in Sec. 3), and define the 
sets A(rk,x) = </>(x, [rfc,oo)). We first observe that L w (x) = f] k A(T k , x). 

PROPOSITION 2 Let x e X Q . If f fl^J /io/cfe ; t/ien 3r = t(x) > such that 
(f>(x,t) £ Y q , Vt > T. 

Proof. We have L u (x)f)Y q = [f] k A(r k , x)] f] Y q = f| fe [^fa, x) f) Now 
let i?(rfe,x) = A(rfc,x)P|y" g . So, for a fixed x £ Xq, the sets £'(rfe,x), k = 
1, 2,3,.. ., is a decreasing sequence of closed subsets of Y q . Since Y q is compact, 
f] fc £'(rfe,x) = implies E(Tk,x) — 0, for all but finitely many fc. Therefore, 
for some if > 1, Vfc > K, E{r k , x) = 0. Hence, 0(x, [i, oo)) f| F g = 0, Vt > t k . 
Therefore, with r = t k , 0(x, f) ^ F 9 , Vt > r. 



We also need the following proposition: 

Proposition 3 Let x e Xq be a point for which there is a t > 0, such that 
0(x, t) ^ y g . TTien t/iere is a S — S x > 0, such that if |x — y| < <5 X and y G Xo, 
t/ien 0(y,r) ^ Y q . 
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Proof. Let |f(x) - f(y)| < C|x - y|. Then |<£(x,r) - 0(y,r)| < |x - y\e CT 
(see, for example, |T^|), so for a fixed r, </>(•, t) is continuous in the first vari- 
able. Let z = 0(x, t). Now since W\l^ is open, there is an e > such that 
B(z,e) C W\Y q . By the continuity of <fi(-, r) at x, there is a (5 > 0, such that 
|0(x, r) — 0(y, t)| = |z — 0(y, t)| < e, whenever |y — x| < S and y e Xo. There- 
fore, c/)(y, t) G B(z, e) C 

From the previous two propositions, we get the following theorem: 

Theorem 2 Let Y q be compact, and f be a continuous function defined on an 
open set W containing Y q satisfying a Lipschitz condition. Further assume 



condition (18) holds for every x G Xq. Then 3r max > independent of x, such 
that for each x G Xq, there is a r(x) with < t(x) < T max and </>(x, r(x)) ^ Y q . 

Proof. Let x G Xq- By Proposition |], there is a r(x) > 0, such that </>(x, t(x)) ^ 
Y" 9 . By Proposition [| there is a i5 x > 0, such that for any y G Xq with 
|y - x| < <S X , <f>(y, r(x)) ^ F g . {S(x, S x ) : x G X } is an open cover of X , 
containing a finite subcover, say, {£>(xi, S Xl ), . . . , £?(x m , S Xm )}. Let r max = 
max{r(xi), . . . , r(x m )}. To check wether this choice of T max is as in the theorem, 
let y G Xq be an arbitrary point. Now y G £> (x^ , 5 Xi ) , for some i with 1 < i < m, 
and by the choice of S Xi , 0(y,r(xi)) g" Y q , concluding the proof. 



Appendix D: Boundedness Results 

In this section, we show that if the initial set F is bounded then the polyhedron, 
P, enclosing -Fjo.A] i as obtained by the algorithm described in Sec. 3 is bounded. 
We assume that F is nonempty. Recall that F is given as the set of points x in 
R" which satisfy the following constraints: 

afx-fe, < 0, i= l,2,...,fc- 1, 
a^x - b k = 0, 

and P is included in the set of points x G W 1 satisfying 

af x - h - k(a£x - b k ) < 0, i = 1, 2, . . . , k - 1, 1 

afcX-6 fc > 0, and \ (20) 

&l*-b k -l k < 0, J 

where (, e R and Zfe > 0. Looking at the constraints, one may visualize the set 
Pi satisfying ( |2p| ) as a prism or a truncated pyramid, with its bottom given by 
( p"9| ) and its top given by 

afx-& i -Z i (a£x-6 fc ) < 0, * = 1, 2, . . . , k - 1, 
a^x -b k -l k = 0. 

We assume that F is nonempty, and wish to show that if F is bounded then 
so is the set P\ , defined as the set of points satisfying (|2(]) . We start with the 
following proposition: 

Proposition 4 Assume that F is nonempty and P\ is not bounded. Then there 
is a point Xo G F and a vector X G 1™, with |A| = 1, such that xo + tX G Pi, for 
all t > 0. 
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Proof. Fix a point xo E F, and suppose Pi is not bounded. So for each 
positive integer k, there is a point Xfc E Pi, such that Ix/j — xq| > k. Let 
A& = fe*~ x "i ■ Now, since F C Pi, Xo € Pi and, by convexity, the entire line 

| x fc x 1 

segment Xo + tXk S Pi, for < t < k. In particular, the points on these line 
segments satisfy (^o|). Now, since for each k, \Xk\ = 1, and the closed unit ball in 
R n is compact, there is a convergent subsequence of A& 's - say, X ki , i = 1, 2, 3, . . . 
- such that Afe i — > A E K n , as i — > oo; since |Afc 4 | = 1, |A| = 1. We show that 
xo + tA € Pi for t E [0, oo). First note that for any positive integer m, and for 
all i > to, xo + iAfei satisfies the constraints (po|), for i 6 [0, fc m ]. Therefore, 
xo + <A satisfies (po|), for t G [0, fc m ], and so Xo + tX € Pi, for < E [0, fe m ]. The 
proposition is concluded by letting to — > oo. 

We now show that such an A, as in the previous proposition, must be parallel 
to the hyperplane passing through P, the normal of which is given by . 

Proposition 5 Ifx E F and a unit vector X E K" are such that x + tX € Pi , 
for all t E [0, oo), then a^A = 0. 

Proof. For a contradiction assume that h — a^A > 0. Now since x + iA satisfies 
(p0|), we must have 

8%(x +tX)-bk-lk<0, 

which holds only if t < < - bk+lk a fc x "^ ; contrary to the hypothesis that, for all 
t E [0, oo), xo + tX E Pi. Similarly, if h = a^A < 0, then the constraint 

af(x + <A)-6 fc >0, 

does not hold for t > 0. Therefore we must have a^A = 0. 

We now show that if F is bounded, then for any unit vector A orthogonal to 
afe, there is an i with 1 < i < k — 1, such that hi = &J X > 0. 

Proposition 6 If F is nonempty and bounded, then for any unit vector X E M. n , 
for which a^A = 0, there is an i, with 1 < i < k — 1, such that &f X > 0. 

Proof. Let xo E F, and suppose that there is an A € K.", with |A| = 1, a^A = 
and &J X < 0, for 1 < i < k — 1. Then for any t > 0, the vector Xo + tX satisfies 
the constraints ([l9|), and so xo + tX E F, contrary to the assumption that F is 
bounded. 

We now combine all the previous propositions to get the following result: 

Theorem 3 If F is nonempty and bounded, then Pi is bounded. 

Proof. If Pi is not bounded, then by Proposition [|, there is an xo E F and a 
unit vector A E W l such that xo + tX E Pi, for any t > 0. By Proposition [|, 
a^A = 0. Now, by the previous proposition, since F is bounded, there is an i, 
with 1 < i < k — 1, such that hi = &f X > 0. But then the constraint 

af (x + tX) -hi- ^(a£(x + tX) -b k )< 0, 



35 



cannot hold for any t with 

af x - bj 

hi 

contrary to the assumption that x + tX e Pi, for any t > 0. Hence P\ is 
bounded. 
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